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Abstract 

In nonlinear voter models the transitions between two states depend in a nonlinear manner 
on the frequencies of these states in the neighborhood. We investigate the role of these 
nonlinearities on the global outcome of the dynamics for a homogeneous network where each 
node is connected to m = 4 neighbors. The paper unfolds in two directions. We first develop 
a general stochastic framework for frequency dependent processes from which we derive 
the macroscopic dynamics for key variables, such as global frequencies and correlations. 
Explicit expressions for both the mean-field limit and the pair approximation are obtained. 
We then apply these equations to determine a phase diagram in the parameter space that 
distinguishes between different dynamic regimes. The pair approximation allows us to identify 
three regimes for nonlinear voter models: (i) complete invasion, (ii) random coexistence, 
and - most interestingly - (iii) correlated coexistence. These findings are contrasted with 
predictions from the mean-field phase diagram and are confirmed by extensive computer 
simulations of the microscopic dynamics. 

PACS: 87.23.Cc Population dynamics and ecological pattern formation, 87.23.Ge Dy- 
namics of social systems 

1 Introduction 

In biological systems, the survival of a species depends on the frequencies of its kin and its foes 
in the environment [3j, |30|]. In some cases, the chance of survival of a certain species improves 



as the frequency of its kind increases, since this might enhance the chance for reproduction or 
other benefits from group interaction. This is denoted as positive frequency dependence. In other 
cases a negative frequency dependence, that is the increase of the survival chance with decreasing 
frequency, is observed. This is the case, when individuals compete for rare ressources. Moreover, 
negative frequency dependence is known to be important for maintaining the genetic diversity 



in natural populations [25l. l39|. 



Frequency dependent dynamics are not only found in biological systems, but also in social and 
economic systems (3, 0, 22, 27, 35, 40, 43, 53]. In democracies, a simple example is a public 



vote, where the winning chances of a party increase with the number of supporters [7j, [12|. In 
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economics, e.g. the acceptance of a new products may increase with the number of its users 
361 ] . In stock markets, on the other hand, positive and negative frequency dependencies may 
interfere. For instance, the desire to buy a certain stock may increase with the orders observed 
from others, a phenomenon known as the herding effect, but it also may decrease, because traders 
fear speculative bubbles. 

In general, many biological and socio-economic processes are governed by the frequency depen- 
dent adoption of a certain behavior or strategy, or simply by frequency dependent reproduction. 
In order to model such dynamics more rigorously (but less concrete), different versions of voter 
models have been investigated. The voter model denotes a simple binary system comprised of N 
voters, each of which can be in one of two states (where state could stand for opinion, attitude, or 
occupation etc.), 9i = {0, 1}. Here, the transition rate w(9\6') from state 0' to state 9 is assumed 
to be proportional to the frequency fg. In this paper, we extend this approach by assuming a 
nonlinear voter model, where the frequency dependence of the transition rate, w{9\9') = n{f) fg, 
includes an additional nonlinearity expressed in terms of the (frequency dependent) prefactor k. 

Linear voter models have been discussed for a long time in mathematics 113. l2lU28i. 12911 . Recently 
they gained more attention in statistical physics 

E 0, H E3, H El, E3, HEJEJH 0, 0, EJ 

because of some remarkable features in their dynamics described in Sect. 12.21 But voter models 
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30, |3J,|3J 



also found the interest of population biologists 19, 23, 

Dependent on how the frequency fg is estimated, one can distinguish global from local voter 
models. In the latter case the transition is governed only by the local frequency of a certain state 
in a given neighborhood. In contrast to global (or mean-field) models, this leads to local effects 
in the dynamics, which are of particular interest in the current paper. If space is represented by a 
two-dimenional lattice and each site is occupied by just one individual, then each species occupies 
an amount of space proportional to its presence in the total population. Local effects such as 
the occupation of a neighborhood by a particular species or the adoption of a given opinion in a 
certain surrounding, can then be observed graphically in terms of domain formation. This way, 
the invasion of species (or opinions) in the environment displays obvious analogies to spatial 
pattern formation in physical systems. 

Physicists have developed different spatial models for such processes. One recent example is 
the so-called "Sznajd model" 0, 0, 0] which is a simple cellular automata (CA) approach to 
consensus formation (i.e. complete invasion) among two opposite opinions (described by spin up 
or down). In 0], we have shown that the Sznajd model can be completely reformulated in terms of 
a linear voter model, where the transition rates towards a given opinion are directly proportional 
to the frequency of the respective opinion of the second-nearest neighbors and independent of 
the nearest neighbors. 

Other spatial models are proposed for game-theoretical interactions among nearest neighbors 
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501 ]. Here, the dynamics are driven by local payoff differences of adjacent players, which 
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basically determine the nonlinearity n(f). Dependent on these payoff differences, we could derive 
a phase diagram with five regimes, each characterized by a distinct spatio-temporal dynamic 



42). The corresponding spatial patterns range from complete invasion to coexistence with large 
domains, coexistence with small clusters, and spatial chaos. 

In this paper, we are interested in the local effects of frequency dependent dynamics in a ho- 
mogeneous network, where each site has m = 4 nearest neighbors. In this case, the nonlinearity 
n{f) can be simply expressed by two constants, ai, «2- This is a special form of a nonlinear 
voter model, which for a\ < 02 < 0.5 also includes majority voting and for a,\ > ai > 0.5 
minority voting. We investigate the dynamics of this model both analytically and by means of 
computer simulations on a two-dimensional stochastic CA (which is a special form of a homoge- 



neous network with m = 4). The latter one was already studied in [30|, in particular there was 
a phase diagram obtained via numerical simulations. In our paper, we go beyond that approach 
by deriving the phase diagram from an analytical approximation, which is then compared with 
our own simulations. 

In Sects. I2l 13.11 we introduce the microscopic model of frequency dependent invasion and demon- 
strate in Sects. 14.11 14.51 the role of ai , «2 by means of characteristic pattern formation. Based 
on the microscopic description, in Sect. I3.2I we derive the dynamics for the global frequency 
x(t), which is a macroscopic key variable. An analytical investigation of these dynamics is made 
possible by pair approximation, Sect. I3.3I which results in a closed-form description for x(t) 
and the spatial correlations C;m(t). In Sect. [5TT1 we verify the performance of our analytical ap- 
proximations by comparing them with averaged CA computer simulations. The outcome of the 
comparison allows us to derive in Sect. I5.2I a phase diagram in the (01,0:2) parameter space, to 
distinguish between two possible dynamic scenarious: (i) complete invasion of one of the species, 
with formation of domains at intermediate time scales, and (ii) random spatial coexistence of 
two species. A third dynamic regime, the nonstationary coexistence of the two species on long 
time scales together with the formation of spatial domains, can be found in a small, but extended 
region that separates the two dynamic regimes mentioned above. We further discuss in Sect. [6] 
that the usual distinctions for the dynamics, such as positive or negative frequency dependence, 
do not necessarily coincide with the different dynamic regimes. Instead, for positive frequency 
dependence, all of the three different dynamic regimes (and the related spatio-temporal patterns) 
are observed. In the Appendix, calculation details for the pair approximation are given. 
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2 Formal Approach to Voter Models 
2.1 Defining the system 

We consider a model of two species labeled by the index a = {0,1}. The total number of 
individuals is constant, so the global frequency x a (or the share of each species in the total 
population) is defined as: 

N = = N + Ni = const. 

a 

X„ = — ; X = Xi = 1 - X (1) 

In the following, the variable x shall refer to the global frequency of species 1. 

The individuals of the two species are identified by the index i G N and can be seen as nodes 
of a network. A discrete value 6{ € {0, 1} indicates whether the node is occupied by species 
or 1. The network topology (specified later) then defines the nearest neighbors ij of node i. 
In this paper, we assume homogeneous networks where all nodes have the same number m of 
nearest neighbors. For further use, we define the local occupation 0^ of the nearest neighborhood 
(without node i) as: 

ii = {0*1)0*21— 1^} ( 2 ) 

A specific realization of this distribution shall be denoted as a, while the function 7/.(er) assigns 
a to a particular neighborhood Of 

£={0"l,cr 2 ,...,<7 TO } 

={0ii=0-l,0i2 =C2,-,0* m = Om) 

For later use, it is convenient to define these distributions also for the nearest neighborhood 
including node i: 

0-i ={0*)0U)0*2)— ,0*'m} = 0* U 

q° ={a, a 1 ,a 2 , ...,<r m } (4) 

jfi (£° ) ={0i = O") 0*i = °"i , 0i 2 = °2 , • • • , 0i m = <7m} 

For m = 4, a denotes a binary string, e.g. {01001}, where the first value a refers to the center 
node and the other values Oj G {0, 1} indicate the particular values of the nearest neighbors. 
The assignment of these values to a particular neighborhood 0° of node i is then described by 

In the voter model described in the following section, the dynamics of 0i is governed by the 
occupation distribution of the local neighborhood, that surrounds each node i. Using a stochastic 
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approach, the probability pi(0i,t) to find node i in state 9i therefore depends in general on the 
local occupation distribution 9 t of the neigborhood (eqn. J2]), in the following manner: 

Pi(M) = S>(0i,&t) (5) 

Hi 

Hence, Pi(9i,t) is defined as the marginal distribution of p(9i, 9 i} t), where 9^ in eqn. (j5J| indicates 
the summation over all possible realizations of the local occupation distribution namely 2 m 
different possibilities. 

For the time dependent change of Pi(9i,t) we assume the following master equation: 



d_ 

dt 



-w{l-9 i \9 i ,tf i )p{9 i ,tf i ,t)\ (6) 

where w(9i\(l — 9i),9 i ) denotes the transition rate for state (1 — 9{) of node i into state 9{ in the 
next time step under the condition that the local occupation distribution is given by 9^. The 
transition rate for the reverse process is w(l — 9i\9i,9 i ). Again, the summation is over all possible 
realizations of 9$, denoted by 9^. It remains to specify the transition rates, which is done in the 
following section. 



2.2 Linear and Nonlinear Voter Models 



Our dynamic assumptions for the change of an individual state 9i are taken from the so-called 
voter model (see also Sect. [T]), abbreviated as VM in the following. The dynamics is given by 
the following update rule: A voter, i.e. a node i £ N of the network, is selected at random and 
adopts the state of a randomly chosen nearest neighbor j. After N such update events, time is 
increased by 1. 

The probability to choose a voter with a given state o from the neighborhood ij of voter i is 
directly proportional to the relative number (or frequency) of voters with that particular state 
in that neigborhood. Let us define the local frequencies in the neighborhood as: 

3=1 

where 5 xy is the Kronecker delta, which is 1 only for x = y and zero otherwise. Then the transition 
rate of voter i to change its state 9i does not explicitly depend on the local distribution 9 it 
but only on the occupation frequency i.e. on the number of nodes occupied by either or 
1 in the neighborhood of size m. Hence, the VM describes a frequency dependent dynamics: 
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the larger the frequency of a given state in the neighborhood, the larger the probability of a 
voter to switch to that particular state if it is not already in that state. I.e. the transition rate 
w(l — 9i\0i = cr, f[), to change state 6 increases only with the local frequency of opposite states, 
fl~ a , in the neighborhood: 

W (i-e t \e l =ajn = 7ft a (8) 

The prefactor 7 determines the time scale of the transitions and is set to 7 = 1. Eqn. j8j describes 
the dynamics of the linear VM because, according to the above update rule, the rate to change 
the state is directly proportional to the frequency. 

The linear (or standard) VM has two remarkable features. First, it is known that, starting from 
a random distribution of states, the system always reaches a completely ordered state, which is 
often referred to as consensus in a social context, or complete invasion in a population biology 
context. As there are individuals with two different states, the complete ordered state can be 
either all or all 1. Which of these two possible attractors of the dynamics is eventually reached, 
depends (in addition to stochastic fluctuations) on the initial global frequency, i.e. x(t = 0). 
It has been shown that, for an ensemble average, the frequency of the outcome of a particular 
consensus state 1 is equal to the initial frequency x(t = 0) of state 1. This second remarkable 
feature is often denoted as conservation of magnetization, where the magnetization is defined as 
M(t) = x\{t) —xo(t) = 2x(t) — 1. Hence, consensus means \M\ = 1. Thus we have the interesting 
situation that, for a single realization, the dynamics of the linear VM is a fluctuation driven 
process that, for finite system sizes, always reaches consensus, whereas on average the outcome 
of the consensus state is distributed as x(0), 

The (only) interesting question for the linear VM is then how long it may take the system 
to reach the consensus state, dependent on the system size N and the network topology. The 
time to to reach consensus, T K) is obtained through an average over many realizations. As the 
investigation of T K is not the focus of our paper (see [2fJ, l2_6|, |49|), we just mention some known 
results for the linear VM : One finds for one-dimensional regular lattices (d = 1) T K oc A 2 and 
for two-dimensional regular lattices (d = 2) T K oc A log A. For d > 2 the system does not always 
reach an ordered state in the thermodynamic limit. In finite systems, however, one finds T K ~ A. 

While the linear VM has some nice theoretical properties, it also has several conceptual disad- 
vantages when applying the model to a social or population biological context. First of all, the 
"voters" do not vote in this model, they are subject to a random (but frequency based) assignment 
of an "opinion", without any choice. Secondly, the state of the voter under consideration does not 
play any role in the dynamics. This can be interpreted in a social context as a (blind) herding 
dynamics, where the individuals just adopt the opinion of the majority. In a population model 
of two competing species, it means that individuals from a minority species may be replaced by 
those from a majority species without any resistance. 
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In order to give voter i at least some weight compared to the influence of its neighbors ij, one can 
simply count its state 9{ into the local frequency ff, i.e. instead of eqs. J2|), J3)) we may use eqn. 
Jl)). Using for voter i the notation 9{ = 9i (i.e j = 0), we can still use eqn. JH]) for the transition 
rates, with the noticable difference that the local frequency ff of eqn. jTJ is now calculated 
from a summation that starts with j = 0. The explicit consideration of 9{ thus has the effect of 
adding some inertia to the dynamics. In fact, extending the summation to j = multiplies the 
transition rate, eqn. ([8]), by a factor m/(m + 1), where m is the number of nearest neighbors. 
I.e., for m = 4 a local configuration a° = {01001} would lead to a transition rate to(l|0) = 0.5 
without the additional inertia, but w(l|0) = 0.4 by counting in the state of voter i. I.e., taking 
into account the state of voter i considerably reduces the transition rate towards the opposite 
state. 

We find it useful for conceptual reasons to include some resistance into the model and therefore 
will use from now on the description which takes the current state of voter i into account. This 
also has the nice advantage that for the case m = 4, which describes e.g. square lattices, we 
avoid stalemate situations, w(l — 9\9) = 0.5. However, we note that the addition of the constant 
resistance does not change the dynamics of the model, as it only adjusts the time scale towards 
a new factor 7' = (m/m + 1) 7. So, keeping m constant and equal for all voters, we can rescale 
7 ' = 1. 

We note that there are of course other ways to give some weight to the opinion of voter i. In 



4a, [47j, we have discussed a modified VM, where voters additionally have an inertia Vi € [0, 1] 



which leads to a decrease of the transition rate to change their state: 

w R (l-9 l \9 l ) = {l-v i )w(l-9 i \9 l ) (9) 

Here w(l — 9\9) is given by the linear VM, eqn. j8]). The individual inertia fj is evolving over 
time by assuming that it increases with the persistence time Tj the voter has been keeping its 
current state. While this intertia may slow down the microscopic dynamics of the VM and thus 
may increase the time to reach consensus, T K , we found the counterintuitive result that under 
certain circumstances a decelerated microdynamics may even accelerate the macrodynamics of 
the VM, thus decreasing T K compared to the linear VM. 

The addition of a nonlinear inertia to the VM, eqn. j9]), is a special case for turning the linear VM 
into a nonlinear one (wheras the fixed resistence would not change the linear VM). In general, 
nonlinear VM can be expressed as 

w(l-9 i \9 i = aj[) = <f) ft* (10) 

where n{f) is a nonlinear, frequency dependent function describing how voter i reacts on the 
occurence of opposite "opinions" in its immediate neighborhood. Fig. [j] shows some possible 
examples which have their specific meaning in a social context. Whereas any function n(f) = 
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const. > describes the linear VM, i.e. a majority voting or herding effect, a decreasing n(f) 
means minority voting, i.e. the voter tends to adopt the opinion of the minority. Nonmonotonous 
n{f) can account for voting against the trend, i.e. the voter adopts an opinion as long as this is 
not already the ruling opinion - a phenomenon which is important e.g. in modeling the adoption 
of fashion. An interpretation of these functions in a population biology context will be given in 
Sect. ED 

In conclusion, introducing the nonlinear response function n(f) will allow us to change the global 
dynamics of the linear VM. Instead of reaching always consensus, i.e. the exclusive domination of 
one "opinion" or species, we may be able to observe some more interesting macroscopic dynamics, 
for example the coexistence of both states. It is one of the aims of this paper to find out, 
under which specifications of k(/) we may in fact obtain a dynamic transition that leads to a a 
structured, but not fully ordered state instead of a completely ordered state. 



A K(f)f 1-0 




against the trend 






y^ linear VM 
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Figure 1: Different nonlinear frequency dependencies for eqn. (|10p . Note the piecewise lin- 
ear functions, as the number of neighbors m and thus the frequencies / have mostly discrete 
values. 

3 Stochastic Dynamics of the Voter Model 
3.1 Microscopic Dynamics 

In order to give a complete picture of the dynamics of the nonlinear VM, we have to derive the 
stochastic dynamics for the whole system of N nodes, whereas eqn. (j6j) gives us "only" the local 
dynamics in the vicinity of a particular voter i. For N voters, the distribution of states is given 

by 

9 = {6 1 ,9 2 ,...,6 N } (11) 

Note that the state space Q of all possible configurations is of the order 2^. In a stochastic 
model, we consider the probability p(@, t) of finding a particular configuration at time t. If t 
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is measured in discrete time steps (generations) and the network is synchronously updated, the 
time-dependent change of p(@, t) is described as follows: 

P (&,t + At) = ^2 P (e,t + At\G',t)p(e',t) (12) 

9' 

where 0' denotes all possible realizations of and p(@,t + At\Q',t) denote the conditional 
probabilities to go from state 0' at time t to at time t + At. Eqn. JT2|) is based on the Markov 
assumption that the dynamics at time t + At may depend only on states at time t. With the 
assumption of small time steps At and the definition of the transition rates 

wmB ', t)=lim E^l±P^l (13) 

v 1 ' ; At-0 At ' 

eqn. ([12]) can be transferred into a time-continuous master equation as follows: 

j t p(@, t) = Y J [w(@\&) p(@', t) - w(G'\G) p(Q, t)] (14) 

G' 

In eqn. (fl4l) . the transition rates depend on the whole distribution 0. However, in the frequency 
dependent dynamics introduced in Sect. 12. ll only the occupation distribution of the local neigh- 
borhood of node i needs to be taken into account. Therefore, it is appropriate to think about 
some reduced description in terms of lower order distributions, such as the local occupation 9 { , 
eqn. J2]). In principle, there are two different ways to solve this task. The first one, the top-down 
approach starts from the global distribution in the whole state space and then uses different 



approaches to factorize p(Q,t). However, a Markov analysis 32J can only be carried out exactly 
for small N, because of the exponential iV- dependence of the state space. Thus, for larger N 
suitable approximations, partly derived from theoretical concepts in computer science need to 
be taken into account. 

In this paper, we follow a second way which is a bottom-up approach based on the local description 
already given in Sect. 12. ll I.e. starting from node i and its local neighborhood, we want to derive 
the dynamics for some appropriate macroscopic variables describing the nonlinear VM. Instead of 
one equation for p(Q, t) in the top-down approach, in the bottom-up approach we now have a set 
of N stochastic equations for pi(8i,t), eqn. (TH]) . which are locally coupled because of overlapping 
neighborhoods, 6 i . In order to solve the dynamics, we need to discuss suitable approximations for 
these local correlations. As we are interested in the macroscopic dynamics, these approximations 
will be done at the macroscopic level. In order to do so, we first derive a macroscopic equation 
from the stochastic eqn. ([6]), which is carried out in the following section. 

3.2 Derivation of the macroscopic dynamics 

The key variable of the macroscopic dynamics is the global frequency x a (t), defined in eqn. JT]). In 
order to compare the averaged computer simulations with results from analytical approximations 
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later in Sect. [5J we first derive an equation for the expectation value (x a ). We do this without 
an explicit determination of the transition rates and wish to emphasize that the formal approach 
presented in Sect. [3] remains valid not just for the voter model, but also for other dynamic 
processes which depend on neighbor interactions (not only nearest neighbors) in various network 
topologies. 

For the derivation of the expectation value we start from the stochastic description given in Sect. 
3.11 where p(@,t) denoted the probability to find a particular distribution 0, eqn. ([!!) ■ at time 
t and O' denoted all possible realizations of eqn. (fill) . On one hand: 

(x a (t)) = 




6' 



and on the other hand: 



1 N 

i=\ 0' 



1 N 

J2pM = °,t) (16) 



N 

i=l 

By differentiating eqn. (fTBT) with respect to time and inserting the master eqn. j6]), we find the 
following macroscopic dynamics for the network: 

N 



i=i 

-w(l-<7|<7,^)p(fl i = <r,^,t)] (17) 

For the further treatment of eqn. (TTT|) ; we consider a specific distribution of states onm+1 
nodes defined by a°. This distribution is assigned to a particular neighborhood of node i by 
V®(<L°) (eqn. Q). Since we are interested in how many times a special realization of a specific 
distribution <7° is present in the population, we define an indicator function 

x(rf l )^xArf i (a°)) = 5 e _ iS _ (18) 

that is 1 if the neighborhood of node i has the distribution a°, and otherwise. Therefore, we 
write the frequency of the n-tuplet a in the population as: 

1 N 

i=l 
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The expectation value is 

(x At)) = ^2x^(0') p(@',t) (20) 
G' 

Inserting eqn. (JX9J) into eqn. (|20l) . we verify that 



i=i e' 



,.i e- (21) 



1=1 



because of the definition of the marginal distribution. Using the identity p(r).,i) = p(a,rj.,t), 
we may rewrite eqn. (|17|) by means of eqn. (|2ip to derive the macroscopic dynamics in the final 
form: 



s! (22) 

a' denotes the 2 m possible configurations of a specific occupation distribution a, eqn. (J3J) ) . In 
the following, we use (x) = (x\) = 1 — {xq). Then, the dynamic for (x) reads: 

±(x(t)) = ^ma') (x 0t At)) 

(23) 

-w(0\l,a') {x lt At)) 

The solution of eqn. (j23j) would require the computation of the averaged global frequencies (x\A 
and (xqA for all possible occupation patterns a, which would be a tremendous effort. Therefore, 
in the next section we will introduce two analytical approximations to solve this problem. In Sect. 
15.11 we will further show by means of computer simulations that these approximations are able 
to describe the averaged dynamics of the nonlinear VM. 



3.3 Mean-Field Limit and Pair Approximation 

As a first approximation of eqn. (|23l) . we investigate the mean-field limit. Here the state of each 
node does not depend on the occupation distribution of its neighbors, but on m randomly chosen 
nodes. In this case the occupation distribution factorizes: 

rn 

(x g o) = (x <r ) l[(x aj ) (24) 
j'=i 
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For the macroscopic dynamics, eqn. (|23|) . we find: 

, m 

- (x(t)) = £ [w(l\0,a>) (1 - <*)) J] <^> 

J =1 . (25) 

-w(p\i,e!) (x)H(x a .) 

3=1 

For the calculation of the (x) we have to look at each possible occupation pattern a for a 
neighborhood m. This will be done in detail in Sect. 14.11 Before, we discuss another analytical 
approximation which solves the macroscopic eqn. (|23p with respect to correlations. This is the 
so-called pair approximation, where one is not interested in the occupation distribution of a whole 
neighborhood a , eqn. (Ill) ) but only in pairs of nearest neigbor nodes, a, a' with a' G {0,1}. 
That means the local neighborhood of nearest neighbors is decomposed into pairs, i.e. blocks of 
size 2 that are called doublets. 

Similar to eqn. (fl9l) . the global frequency of doublets is defined as: 

. jv m - 

V^EE-xl^Mi,^') (26) 
i=i j=i 

The expected value of the doublet frequency is then given by {x a ,a') m the same way as in eqn. 
(|2"0"1) , We now define the correlation term as: 

c ° w> - ~wr (27) 

neglecting higher order correlations. Thus c a \ a i can be seen as an approximation of the conditional 
probability that a randomly chosen nearest neighbor of a node in state a' is in state a. Using 
the above definitions, we have the following relations: 

(av) c a \ a , = (x a ) <V|o. ; 52 C °'W = 1 ( 28 ) 

a'e{0,l} 

For the case of two species a € {0,1}, C\\\ and c |o are the inter-species correlations, while 
c 1 |g and c n denote the intra-species correlations. Using (x) = (xi), these correlations can be 
expressed in terms of only c^h and (x) as follows: 

Co|l =1 - c l|l 

C H0 " i _ ( X ) (29) 
_l-2(x) + (x) ci|i 

c °|° - 
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Now, the objective is to express the global frequency of a specific occupation pattern (x a o), 
eqn. (|20|) . in terms of the correlation terms c a \ a i. In pair approximation, it is assumed that the 
states 9i j are correlated only through the state 9i and uncorrelated otherwise. Then the global 
frequency terms in eqn. (|22j) can be approximated as follows: 

rn 

(v} = W[[si' ( 30 ) 

i=i 

For the macroscopic dynamics, eqn. ([23]) . we find in pair approximation: 

, m 

Jt <*(*)) = Y, [<"(i|0.2')(i - M)IK-k 

(3D 

-w(0|1,ct') (x) n c ^l(i-^) 
i=i 

Note that the c a .\ a can be expressed in terms of c^i by means of eqn. (|29j) . Thus, eqn. (|3lT) now 
depends on only two variables, (x) and em. In order to derive a closed form description, we need 
an additional equation for em - That can be obtained from eqn. (f2Tj) : 

^r = -V)dt {x) + V)dt {xi > l) (32) 

Eqn. J321) also requires the time derivative of the global doublet frequency Even in their 

lengthy form, the three equations for (x), cm, (x^i) can easily be solved numerically. This gives 
the approach some computational advantage compared to averaging over a number of microscopic 
computer simulations for all possible parameter sets. 

Although the approach derived so far is quite general in that it can be applied to different network 
topologies and neighborhood sizes, specific expressions for these three equations of course depend 
on these. Therefore, in the Appendix, these three equations are explicitly derived for a 2d regular 
lattice with neighborhood m = 4 using the specific transition rates introduced in the next section. 
In Sect. 15.11 we further show that the pair approximation yields some characteristic quantities 
such as (x(t)) for the 2d regular lattice in very good agreement with the results of computer 
simulations. 

4 Invasion versus Coexistence 
4.1 Nonlinear Response Functions 

So far, we have developed a stochastic framework for (but not restricted to) nonlinear voter 
models in a general way, without specifying two of the most important features, namely (i) the 
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network topology which defines the neighborhood, and (ii) the nonlinearity n{f) which defines 
the response to the local frequencies of the two different states. For (i), let us choose a regular 
network with m = 4, i.e. each voter has 4 different neighbors. We note explicitly that our 
modeling framework and the general results derived hold for all homogeneous networks, but for 
the visualization of the results it will be most convenient to choose a regular square lattice, where 
the neighbors appear next to a node. This allows us to observe the formation of macroscopic 
ordered states in a more convenient way, without restricting the general case. Eventually, to 
illustrate the dynamics let us now assume a population biology context, where each node is 
occupied by an individual of either species or 1. The spreading of one particular state is then 
interpreted as the invasion of that respective species and the local disappearence of the other one, 
while the emergence of a complete ordered state is seen as the complete invasion or domination 
of one species together with the extinction of the other one. 

Keeping in mind that we also consider the state of node i itself, we can write the possible 
transition rates, eqn. (flu]) , for the neighborhood ofra = m+ l = 5 and 9{ = a in the following 
explicit way (cf also [30] ) : 



J I 


Al-a) 
J i 






1 







a 


4/5 


1/5 




«i 


3/5 


2/5 




a -2 


2/5 


3/5 




= 1 — 0*2 


1/5 


4/5 


«4 


= 1 — Oi\ 





5/5 


Or, 


= l-a 



Eqn. (|33p means that a particular node i currently in state 9{ = a, or occupied by an individual 
of species a where a is either or 1, will be occupied by an individual of species (1 — a) with 
a rate w(l — 9i\9i = a, f[) that changes with the local frequency ff in a nonlinear manner. The 
different values of a n denote the products ^{f)fl~ a for the specific values of / given. I.e., the 
a n define the piecewise linear functions shown in Fig. [TJ 

The general case of six independent transition rates a n (n = 0,...,5) in eqn. ([33]) can be reduced 
to three transition rates an, «i, «2 by assuming a symmetry of the invasion dynamics of the two 
species, i.e. a% + 013 = 1, ot\ + 04 = 1 and an + as = 1. Further, assuming a pure frequency 
dependent process, we have to consequently choose an = 0, because in a complete homogeneous 
neighborhood, there is no incentive to change to another state (there are no other species around 
to invade). 

We recall that if the transition rates a±, ai are directly proportional to f^ a \ i.e. a\ = 0.2 
and 0C2 = 0.4, this recovers the linear VM, eqn. (jHJ). (Note that without the resistence of node 
i discussed in Sect. 12.21 the linear voter point would read as a% = 0.25 and 02 = 0.5 instead.) 
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Dependent on the relation of the two essential parameters ai, a 2 , we also find different versions 
of nonlinear VM, which have their specific meaning in a population biology context: 

(pf) < ax < a 2 < (1 - a 2 ) < (1 - ai) < 1 
(nf) 1 > ai > a 2 > (1 - a 2 ) > (1 - «i) > 
(pa) < qi < a 2 , 012 > (1 - a 2 ), 

(1 -a 2 ) < (1 - ai) < 1 
(na) 1 > ai > a 2 , a 2 < (1 - a 2 ), 

(1 -a 2 ) > (1 -ai) > (34) 

Note, that the parameters ai,...,04 can be ordered in 24 different ways. These reduce to 8 
inequalities under the conditions 03 = 1 — a 2 and «4 = 1 — ai. In eqn. (I3"4"l) . (pf) means 
(pure) positive freqency dependent invasion, where the transition rate increases with an increasing 
number of individuals of the opposite species (1 — a) in the neighborhood, and (nf) means (pure) 
negative freqency dependent invasion because the transition rate decreases. The two other cases 
describe positive (pa) and negative (na) allee effects [30J]. These regions are described by 3 
inequalities each, all of which show the same relative change in parameter values, if going from 
a± to 0:4. Similar to the drawings in Fig.[T]this can be roughly visualized as an up-down- up change 
in region (pa) and a down-up-down change in region (na). The different parameter regions are 
shown in Fig. [2j On a first glimpse, one would expect that the dynamics as well as the evolution 
of global variables may be different in these regions. Thus, one of the aims of this paper is to 
investigate whether or to what extent this would be the case. 

1.0 




0.2 0.4 0.6 0.8 1.0 



Figure 2: Four different parameter regions for frequency dependent invasion, according to eqn. 
(|34l) . The linear voter point is indicated by • 
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4.2 Mean-Field Analysis 

In order to find out about the influence of the nonlinear response function «(/), which is specified 
here in terms of a\, a%, let us start with the mean- field approach that lead to eqn. (|25l) . As we 
outlined in Sect. 13.31 the calculation of the (x aj ) in eqn. ([25]l requires to look at each possible 
occupation pattern a for a neighborhood m, for instance, a = {0010}. The mean- field approach 
assumes that the occurence of each 1 or in the pattern can be described by the global frequencies 
x and (1 — x), respectively (for simplicity, the abbreviation x = (x) will be used in the following). 
For the example of string a = {0010} we find F] (x aj ) = x(l — x) 3 . The same result yields 
for a = {0100} and for any other string that contains the same number of 1 and 0, i.e. there 
are exactly ( 4 ) different possibilities. For strings with two nodes of each species, ( 2 ) times the 
contribution n( x o-,) = x 2 (l — x) 2 results, etc. Inserting eqn. ([33]) for the transition rates, we 
find with an = the equation for the mean- field dynamics: 



-y- = (1 - x) 4aix(l - x) 3 + 6a 2 x 2 (l - x) 2 



dt 



+ 4(1 - a 2 )x 3 (l - x) + (1 - ai)x 4 

1 - ai)(l - x) 4 + 4(1 - ai)x(l - x) 3 
+ 6a 2 x 2 (l - x) 2 + 4aix 3 (l - x) 



(35) 



The fixed points of the mean-field dynamics can be calculated from eqn. (|35"1) using x = 0. We 
find: 

X W = ; x (2) = 1 ; x (3) = 0.5 

(4,5) 



x( 4 ' 5 ) = 0.5±VAM (36) 
/3\ = a 2 + 1.5«i — 0.7; = «2 — 0.5ai — 0.3 

The first three stationary solutions denote either a complete invasion of one species or an equal 
share of both of them. "Nontrivial" solutions, i.e. a coexistence of both species with different 
shares of the total population, can only result from x^ 4,5 ), provided that the solutions are (i) 
real and (ii) in the interval {0, 1}. The first requirement means that the two functions /3 2 are 
either both positive or both negative. The second requirement additionally results in a\ < 0.2 if 
o~2 > 0.4 and a\ > 0.2 if a 2 < 0.4. This leads to the phase diagram of the mean-field case shown 
in Fig. [3l 

In order to verify the stability of the solutions, we have further done a perturbation analysis (see 
also Sect. I4.4I) . The results can be summarized as follows: 

• In the regions a and b of the mean-field phase diagram, Fig. OH x stat = and x stat = 1 are 
the only stable fixed points of the dynamics, while x = 0.5 is an unstable fixed point (cf. 
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Figure 3: Phase diagram of the invasion dynamics in the mean- field case. For the different areas 
see text. The functions f3\ = and = are given by eqn. jM]). The areas in lighter gray 
indicate imaginary solutions of x^ 4,5 \ eqn. (|36l) while the areas in darker gray (color: yellow) 
indicate solutions of x( 4 ' 5 ) outside the {0, 1} interval. The linear voter point is indicated by •. 
Note that the physically relevant solutions x 1,2 ' 3 are the same in the lighther gray and darker 
gray (color: yellow) areas, however their stability is different in (a,b) and (d,e). 




Figure 4: Bifurcation diagram of the stationary solutions dependent on a\ and ct2- (top) a.\ = 0.1, 
(bottom) «2 = 0.1. The solid lines refer to stable solutions, the dashed lines to unstable ones. 
The notations a-f refer to the respective areas in the phase diagram, Fig. [3j 

also Fig. H]top). Species 1 with x(t = 0) < 0.5 will most likely become extinct, while it will 
remain as the only survivor for x(t = 0) > 0.5. Thus, the region (a, b) can be characterized 
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as the region of invasion. 

• In region c, the mean-field limit predicts the three stable fixed points 0, 1 and 0.5. The 
attractor basin for 0.5 is the largest as Fig. IDJtop) indicates. The separatices are given 
by the unstable solutions a^ 4 ' 5 -*, eqn. (|36l) . In this parameter region, the mean-field limit 
predicts either coexistence of both species with equal shares, or invasion of one species, 
dependent on the initial condition x(t = 0). 

• In the regions d and e, only one stable fixed point x stat = 0.5 can be found, while the 
solutions and 1 are unstable (cf. also Fig. Hbottom). Thus, the mean-field approach 
predicts the coexistence of both species with equal share. 

• Finally, in region / the solutions 0, 1 and 0.5 are unstable fixed points , but the two 
remaining solutions x^ 4,5 ), eqn. (|36|) are stable fixed points (cf. Fig. HJaottom). Thus, this 
region is the most interesting one, since it seems to enable "nontrivial" solutions, i.e. an 
asymmetric coexistence of both species with different shares. We note again, that this is 
a prediction of the mean-field analysis. At the intersection of regions / and a, these two 
solutions approach and 1, while at the intersection of regions / and e they both converge 
to 0.5. 

We will compare these mean-field predictions both with computer simulations and analytical 
results from the pair approximations later in this paper. Before, in Sects. I4~3l 14.41 we would like 
to point to some interesting (01,02) combinations in this phase diagram where the mean-field 
analysis does not give a clear picture of the dynamics. 

4.3 Deterministic limit 

The first set of interesting points are (01,02) combinations of values and 1, such as (01,02) = 
(0, 0) etc. These cases are special in the sense that they describe the deterministic limit of the 
nonlinear voter dynamics. Whereas for < o n < 1 always a finite probability exist to change 
to the opposite state, for (0, 0) the state of node i never changes as long as at least half of 
the nearest neighbor nodes are occupied by the same species. On the other hand, it will always 
change if more then half of the neighboring nodes are occupied by the other species. This refers to 
a deterministic positive frequency invasion process. Similarly, a deterministic negative frequency 
invasion process is described by (1, 1). 

The deterministic dynamics, as we know from various other examples, may lead to a completely 
different outcome as the stochastic counterpart. In order to verify that we have conducted com- 
puter simulations using a cellular automaton (CA), i.e., a two-dimensional regular lattice with 
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periodic boundary conditions and synchronous update of the nodes. The latter one can be ar- 
gued, but we verified that there are no changes in the results of the computer simulations if the 
sequential update is used. The time scale for the synchronous update is defined by the number 
of simulation steps. If not stated otherwise, the initial configuration is taken to be a random dis- 
tribution (within reasonable limits) of both species, i.e. initially each node is randomly assigned 
one of the possible states, {0, 1}. Thus, the initial global frequency is x(t = 0) = 0.5. Fig. [5] shows 
snapshots of computer simulations of the deterministic dynamics taken in the (quasi-)stationary 
dynamic regime. 




(c) (0,0) (d) (1,0) 



Figure 5: Spatial snapshots of a deterministic frequency dependent invasion process for different 
combinations of (01,02). The pictures are taken after t = 10 2 time steps, except for (a), where 
t = 10 4 . Lattice size: 80 x 80. In all pictures black dots refer to species 1. 



If we compare the snapshots of the deterministic voter dynamics with the mean-field prediction, 
the following observations can be made: 



1. A spatial coexistence of both species is observed for the (01,02) values (0,0), Fig. 5(c 
(1,0), Fig. |5(d)[ (1,1), Fig. |5(b)| where the global frequency in the stationary state is 
x stat = 0.5. This contradicts with the mean-field prediction for (0,0), which is part of 
region (a) and thus should display complete invasion. 
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2. A complete 

(ai,a 2 ) = (0,1), Fig. [5(a 



invasion of one species is observed for 
This would agree with the mean-field prediction of ei- 



ther coexistence or invasion. A closer look at the bifurcation diagram, Fig.[H(top), however 
tells us that for the given intial condition x{t = 0) = 0.5 the stationary outcome should 
be coexistence, whereas the deterministic limit shows always invasion as was verified by 
numerous computer simulations with varying initial conditions. 

3. For the deterministic frequency dependent processes (0,0), (1,0) the spatial pattern be- 
comes stationary after a short time. For the negative frequency dependence (1, 1) the 
pattern flips between two different configurations at every time step. So, despite a constant 
global frequency x sta,t = 0.5 in the latter case, local reconfigurations prevent the pattern 
from reaching a completely stationary state, but it may be regarded as (quasi-) stationary, 
i.e. the macroscopic observables do not change but microscopic changes still occur. 

4. In both - positive and negative - deterministic frequency dependence cases, individuals 
of the same species tend to aggregate in space, albeit in different local patterns. For the 

or even 



positive frequency dependence, we see the occurence of small clusters, Fig. 5(c 



complete invasion, Fig. 5(a), based on the local feedback within the same species. For 
the negative frequency dependence however we observe the formation of a meander-like 
structure that is also known from physico-chemical structure formation [4l| . It results from 
the antagonistic effort of each species to avoid individuals of the same kind, when being 
surrounded by a majority of individuals of the opposite species. 



In conclusion, the mean-field analysis given in this section may provide a first indication of how 
the nonlinearities may influence the voter dynamics. This, however, cannot fully extended to the 
limiting cases given by the deterministic dynamics. 



4.4 Perturbation Analysis of the Linear Voter Point 

The other interesting (01,02) combination is the linear voter point (0.2,0.4) where all different 
regions of the mean-field phase diagram, Fig. [3], intersect. Inserting (0.2,0.4) into eqn. (|35l) 
yields dx/dt = regardless of the value of x, i.e. x(t) = x{t = 0) for all initial conditions. This 
important feature of the linear VM was already discussed in Sect. 12.21 We recall that, while 
on the one hand the microscopic realizations always reach consensus (complete invasion of one 
species) in the long term, on the other hand an averaged outcome over many realizations shows 
that the share of the winning species is distributed as x(t = 0). To put it the other way round, 
the mean-field limit discussed failes here because it does not give us any indication of the fact 
that there is a completely ordered state in the linear VM. The averaged outcome, for example 
x = 0.5, can result both from complete invasion of species 1 (50 coexistence of the two species. 
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Both of these outcomes exist in the immediate neighborhood of the linear voter point as Fig. [3] 
shows. In order to get more insight into this, we will later use the pair approximation derived 
in Sect. 13.31 Here, we first follow a perturbation approach, i.e. we add a small perturbation to 
the solution describing the macroscopic ordered state of complete invasion (consensus). In terms 
of the nonlinear response function expressed by the a n in eqn. (1331) . this means a nonzero 

value of a.Q = e, i.e. a small parameter indicating the perturbation. With this, we arrive at a 
modified mean-field equation: 



dt £ 



:i -x) 5 -x 5 



where dx/dt is given by the nonperturbated mean-field eqn. (|35j) and the index p shall indicate 
the presence of the perturbation e. Consequently, this changes both the value of the fixed points, 
previously given by eqn. (f36| and their stability. Instead of a complete analysis in the (e, a±, 02) 
parameter state, we restrict the investigations to the vicinity of the linear voter point (ai, 02) = 
(0.2,0.4) where dx/dt = 0. Eqn. (|3*7)l then returns only one real stationary solution, Xp tat = 0.5, 
which is independent of e and stable. Consequently, any finite perturbation will destroy the 
characteristic feature of reaching an ordered state in the linear VM, i.e. complete invasion, and 
leaves only coexistence of both species as a possible outcome. This is little surpising because 
adding an an > to the dynamics transforms the former attractor x — > 0, 1 into a repellor, i.e. it 
prevents reaching the ordered state. More interesting the question is, how the perturbated linear 
voter dynamics looks in detail. This is investigated in the next section by means of computer 
simulation, and in Sect. I5.ll by means of the pair approximation approach. 



4.5 Computer simulations of the perturbated CA 

For further insight into the dynamics of the nonlinear VM, we perform some computer simulations 
using the CA approach already described in Sect. 14. 3( Is important to notice that we have chosen 
different sets of the parameters ai, «2 from the region of positive frequency dependence, as defined 
in Fig. [21 I.e., the transition towards the opposite state strictly increases with the number of 
neighbors in that state (majority voting). So one would naively expect a similar macroscopic 



dynamics in that region as done in [30(| This however is not the case as the following simulations 
indicate. A thorough analysis is presented in Sect. 15.21 

In order to study the stability of the global dynamics for the different a\ , 0,2 settings in the vicinity 
of the linear voter point, we have added a small perturbation an = £• As the investigations of 
Sect. 14.41 have indicated, we should no longer expect consensus for the perturbated linear VM, 
but some sort of coexistence. In fact, we observe an interesting nonstationary pattern formation 
we call correlated coexistence. Fig. [6] (obtained for another range of parameters a%, 02) shows 
an example of this. We find a long-term coexistence of both species, which is accompanied by 
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a spatial structure formation. Here, the spatial pattern remains always nonstationary and the 
global frequency randomly fluctuates around a mean value of x = 0.5, as shown by Fig. [71 




Figure 6: Spatial snapshots of a positive frequency dependent invasion process with a\ = 0.24, 
«2 = 0.30, e = 10~ 4 (nonstationary, correlated coexistence), (a) t = 10 1 , (b) t = 10 2 (c) t = 10 3 , 
(d) t = 10 4 . Note that a simulation using the parameters of the linear VM, a\ = 0.2, «2 = 0.4, 
would look statistically similar, and also the fluctuations of the global frequencies shown in Fig. 
[TJare quite similar. 



In more specific terms, the regime defined as 'correlated coexistence' is a paramagnetic phase 
with finite domain length, typical of partially phase-separated systems. We mention that such 
a regime was also observed in some related investigations of the VM and other nonlinear spin 
models with Ising behavior [la, [13]. Also, a similar transition was observed in the Abrams- 



Strogatz model [1], where the transition rate is a power a of the local field. The stability of the 
solutions then changes at a = 1, from coexistence for a < 1 to dominance for a > 1. 

In order to find out about the range of parameters in the nonlinear VM resulting in the quite 
interesting phenomenon of correlated coexistence, we have varied the parameters a%, ai within 
the region of positive frequency dependence. Fig. factually shows results from a set picked from 
region (f) in Fig. EJ where the mean-field analysis predicts an asymmetric coexistence of both 
species. Obviously, the nonstationarity results from the perturbation e. 

However, for other sets a±, ai in the positive frequency dependence region the perturbation does 
not prevent the system from reaching a global ordered state, i.e. invasion of one species as Fig. [8] 
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Figure 7: Global frequency of species 1 vs. time (a) for the linear voter model, ot\ = 0.2, «2 = 0.4, 
e = 10" 4 , and (b) for the same setup and parameters as in Fig. [6] The initial frequency is 
x(t=0)=0.5 for both runs. 

verifies. This process is accompanied by a clustering process and eventually a segregation of both 
species indicated by the formation of spatial domains. Fig. [9] depicts the evolution of the global 
frequency x(t) of species 1 for different intitial frequencies x(t=Q). In every case, one species 
becomes extinct. For x(t=Q) > 0.5 species 1 is the most likely survivor, while for x(t = 0) < 0.5 
it is most likely to become extinct. For x(t = 0) = 0.5, random events during the early stage 
decide about its outcome. 

On the other hand, the perturbation also does not induce an ordered state as the random co- 
existence in Fig. [10] shows, which was again obtained from parameter settings in the region of 
positive frequency dependence. So, we conclude that computer simulations of positive frequency 
dependent processes show three different dynamic regimes (dependent on the parameters ai, 
02): (i) complete invasion, (ii) random coexistence, and (iii) correlated coexistence. While for (i) 
and (ii) the outcome is in line with the mean-field prediction shown in Fig. HJ this does not im- 
mediately follows for (iii). So, we are left with the question whether the interesting phenomenon 
of correlated coexistence is just because of the perturbation of some ordered state, or whether it 
may also exist inspite of e. 

We just add that for negative frequency dependent invasion, eqn. J34J, the the spatial pattern 
remains random, similar to Fig. [10] Furthermore, regardless of the initial frequency x(t = 0) , on a 
very short time scales, a global frequency x stat =0.5 is always reached. That means we always find 
coexistence between both species. We conclude that for negative frequency dependence x stat =0.5 
is the only stable value, which is in agreement with the mean-field prediction, whereas for positive 
frequency dependence the situation is not as clear. 
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(d) 



Figure 8: Spatial snapshots of a positive frequency dependent invasion process with a\ 
a 2 = 0.3, e = 10" 4 (complete invasion), (a) t = 10 1 , (b) t = 10 2 (c) t = 10 3 , (d) t = 10 4 . 



0.1, 




500 1000 1500 2000 2500 
t 



Figure 9: Global frequency of species 1 vs. time for the same setup and parameters as in Fig. 
The initial frequencies x(t=0) of the four different runs are: (a) 0.6, (b) 0.5, (c) 0.5, (d) 0.4. 



5 Derivation of a phase diagram 



To answer the question what ranges of a\, ai eventually lead to what kind of macroscopic 
dynamics, we now make use of the pair approximation already derived in Sect. 13.31 as a first 
correction to the mean-field limit. Here, we follow a two-step strategy: First, we investigate how 



24/19] 



Frank Schweitzer, Laxmidhar Behera: 
Nonlinear Voter Models: The Transition from Invasion to Coexistence 

European Physical Journal B, vol 67 (2009). 



Figure 10: Spatial snapshots of a positive frequency dependent invasion process with ol\ = 0.3, 
ct2 = 0.4, e = 10~ 4 (random coexistence). The snapshot shown at t = 10 5 is statistically 
equivalent to the initial random state. 

well the pair approximation, eqs. j38j), (l43j) . (|44"1) of the macroscopic dynamics, eqn. (|23j) . predict 
the global quantities (x) and cih. In order to specify the network topology, we use again the CA 
described above. Second, we use the pair approximation to derive a phase diagram in the case of 
local interaction. Eventually, we test whether these findings remain stable against perturbations 
of the ordered state. All predictions are tested by comparison with computer simulations of the 
microscopic model, from which we calculate the quantities of interest and average them over 50 
runs. 

5.1 Global frequencies and spatial correlations 

Here we have to distinguish between the three different dynamic regimes already indicated in 
Sect.l4~5l 

Regime (i), complete invasion, is characterized by fixed points of the macroscopic dynamics of 
either {{x) ,cin} = {1, 1} or {(x) , ci^} = {0,0}. The CA simulations as well as also the pair 
approximation of the dynamics quickly converge to one of these attractors, dependent on the 
initial conditions. 

Regime (ii), random coexistence, has only one fixed point, {(x) , cui} = {0.5,0.5}, to which the 
CA simulations quickly converge. The pair approximation converges to (x) = 0.5 after some initial 
deviations from the CA simulation, i.e. it relaxes on a different time scale (t > 40), but is correct 
in the long run. The approximation of the local correlation cm shows some deviations from the 
predicted value, cm = 0.5. We have tested the case of random coexistence for various parameter 
values and found values for cm between 0.4 and 0.6. The discrepancy is understandable, since in 
the case of long-term coexistence some of the spatial patterns flip between two different random 
configurations with high frequency. Thus, while the global frequency settles down to 0.5, the 
microscopic dynamics is still nonstationary. 
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Regime (iii), correlated coexistence, the most interesting one, is chacterized by an average global 
frequency of (x) = 0.5 again, however the existing local correlations lead to a much higher 
value of C]h > 0.6. This is shown in Fig. [TTJ where we find cm ~ 0.7 from the CA simulations 
and C!|! « 0.65 from the pair approximation. I.e., for the case of spatial domain formation the 
long-range correlations can be well captured by the pair approximation, whereas this was less 
satisfactory for the short-range correlations of the random patterns. 

i | i | i | i | i I I i | i | i | i | i 




4000 6000 8000 10000 



t t 



Figure 11: Global frequency (x) with min-max deviations (top) and spatial correlation cm 
(bottom) for the case of correlated coexistence. The two curves shown in the bottom part result 
from averaging over 50 CA simulations (black dotted line) and from the pair approximation 
(green solid line) The parameters are as in Figs. [6j El (top) 

5.2 Determining the phase boundary 

The insights into the dynamics of the nonlinear VM derived in this paper are now summarized in 
a phase diagram that identifies the different parameter regions (01,02) f° r the possible dynamic 
regimes identified in the previous sections. In order to find the boundaries between these different 
regimes, we carried out CA simulations of the complete parameter space, < (01,02) < 1. 
Precisely, for every single run the long-term stationary values of x and cm were obtained and 
then averaged over 50 simulations. As described above, the three different regimes could be clearly 
separated by their {(x) ,cm} values, which were used to identify the phase boundary between 
the different regimes. The outcome of the CA simulations is shown in the phase diagram of Fig. 
ITU (top) and should be compared with Fig. El which results from the mean-field analysis in Sect. 
13.31 and thus neglects any kind of correlation. 
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Figure 12: Phase diagram of the nonlinear voter model: (top) CA simulations, averaged over 
50 runs (cm = 0.7), (bottom) pair approximation. Phase boundaries for e = resulting from 
cui = 0.65 (upper limit): left - red solid line, right - black solid line, for comparison phase 
boundaries resulting from cm = 0.60 (lower limit): left - red solid line (identical with c\u = 0.65), 
right - red solid line on the far right. Dashed red lines indicate the shift of the phase boundaries 
for C]h = 0.60 if e = 0.02. The linear voter point (0.2,0.4) is indicated by •. Further, O marks 
those three parameter sets from the positive frequency dependence region where CA simulations 
are shown in Figs. I8l l6l [TUI (see also global frequency (x) and local correlation am in Figs. [TJ EE 
ITTT) . The straight dashed lines mark the different parameter areas given in eqn. ([Ml) which are 
also shown in Fig. [2] 

Instead of the six regions distinguished in Fig. [3] in the local case we can distinguish two different 
regions divided by one separatrix: The parameter region left of the separatrix refers to the 
complete invasion of one of the species, with high local correlations during the evolution. In 
the CA simulations, we observe the formation of domains that grow in the course of time until 
exclusive domination prevails. Asymptotically, a stationary pattern is observed, with (x) = 1 (or 
0) and Cx|i = 1 (or 0) (cf the simulation results show in Figs. [HI EI) • I.e., the system converges 
into a 'frozen' state with no dynamics at all. The region to the right of the separatrix refers to 
random coexistence of both species with (x) = 0.5 and no local correlations, i.e. Ci|i = 0.5. In 
the CA simulations, we observe nonstationary random patterns that change with high frequency 
(cf. Fig. [TO]). 

Both regions are divided by a separatrix. As shown in Fig. [12] (top), the separatrix is divided into 
two pieces by the linear voter point, (0.2,0.4). Above that point, the separatrix is very narrow, 
but below the voter point it has in fact a certain extension in the parameter space. Looking at 
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the dynamics on the separatrix, we find that the mean frequency is (x) = 0.5 both above and 
below the voter point (see also Fig. [TTJ top). The local correlation cm = 0.7 holds within the 
whole extended area of the separatrix which identifies it as the region of correlated coexistence 
(see also Fig. [Til bottom). 

The question is how well this phase diagram can be predicted by using the pair approximation 
of the dynamics, described by eqs. ([38]) . (pHl) . For a comparison, the coupled equations were 
numerically solved to get the asymptotic solutions for {(x) ,cm} (which looks complicated but 
is numerically very fast and efficient). In order to distinguish between the three regimes, we 
have to define a critical cm value for the correlations. Whereas in the CA simulations cm = 0.7 
indicated a correlated nonstationary coexistence, this value was never reached using the pair 
approximation (see Fig. [TTJ bottom), so cm < 0.65 could be regarded as an upper limit for the 
case of correlated coexistence. Random coexistence, on the other end, yielded cm = 0.5 in the 
CA simulations and a value between 0.4 and 0.6 in the pair approximation. So, given suitable 
initial conditions, we can regard cm > 0.6 as a lower limit, to identify correlated coexistence. 
The results are shown in Fig. [12] (bottom) which shall be compared with the phase diagram 
above (Fig. CEJ top). 

Fig. [12] (bottom) shows the influence of the cm threshold value. With the lower limit, we find a 
quite broad region of correlated coexistence, which for example also includes the point (0.3, 0.4) 
for which a random coexistence in the CA simulations was shown in Fig. [10] Using the upper 
limit cm = 0.65 results in a much smaller region of correlated coexistence. Comparing this 
with the CA simulations above, we can verify that the pair approximation correctly predicts the 
extended region below the voter point and also shows how it becomes more narrow above the 
voter point. One should note that the left border of the separatrics is not affected by the threshold 
value, whereas the right border shifts considerably. The left border also contains the voter point 
(independent of the cm threshold value), for which a complete invasion can be observed. 

Therefore, it is quite interesting to look into changes of the phase diagram if additional pertur- 
bations are considered (see also Sect. 14.4(1 . Fig. [12] (bottom) shows (for the threshold cm = 0.6) 
that this does not affects the existence of the three dynamic regimes and most notably of the 
extended separatrix below the voter point, but only shifts the boundaries toward the left, depen- 
dent on the value of an = £ (this can be also verified for cm = 0.65 but is omitted here, to keep 
the figure readable). Thus, the consideration of perturbations in the phase diagram reveals that 
it is indeed the nonlinearity in the voter model which allows for the interesting phenomenon of 
the correlated coexistence and not just the perturbation. 

A closer look into Fig. [12] (bottom) also shows that in the perturbated phase diagram the voter 
point no longer lies on the boundary towards the region of complete invasion but clearly within 
the region of correlated coexististence. This is in agreement with the findings in Sect. 14.41 which 
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showed that for the linear VM complete invasion is an unstable phenomenon and changes into 
correlated coexistence for finite e. 



6 Discussion and Conclusions 

In this paper, we investigated a local model of frequency dependent processes, which for example 
models the dynamics of two species {0, 1} in a spatial environment. Individuals of these species 
(also called 'voters') are seen as nodes of a network assumed as homogeneous in this paper (i.e. 
all nodes have the same number of neighbors, m). The basic assumption for the microscopic 
dynamics is that the probability to occupy a given node with either species or 1 depends on 
the frequency of this species in the immediate neighborhood. Different from other investigations, 
we have counted in the state of the center node as well (see Sect. 12.21) and have further considered 
a nonlinear response of the voters to the local frequencies. 

Studies of a nonlinear version of the traditional voter model (without counting the state of the 



central node and with sequential dynamics) have already been analyzed before. Ref. [30||, as 
pointed out before, is closest to our investigations, but restricted itself to the mean-field analysis 
and computer simulations of the 2d case, to obtain a phase diagram similar to Fig. [12] (top). 



321 ] . on the other hand, have provided a Markov analysis which is restricted only to very small 



CA. The two-parameter model in (15J] is for y = 1 a nonlinear voter model that exhibits at the 
voter point (x = 1/2) a transition from a ferromagnetic phase, i.e. invasion, for x > 1/2 to a 
paramagnetic phase (correlated coexistence) for x < 1/2. Also the case of the 'perturbated linear 
voter model' is included in the model, for x = 1/2 and y < 1. Similar results are also presented in 



181 ] . which points out relations to random branching processes, and in [111 ], where the emphasis 



is on investigations of the interface density, to describe the coarsening process. A recent paper 



521 ] also shows for spin systems with two symmetric absorbing states (such as the VM) that the 
macroscopic dynamics only depends on the first derivatives of the spin-flip probabilities. 

In our paper, we set out for a formal approach that allows to derive the dynamics on different 
levels: (i) a stochastic dynamics on the microlevel, which is used for reference computer simula- 
tions but also allows a derivation of the (ii) macroscopic dynamics for the key variables, given in 
terms of differential equations. This macroscopic dynamics is then analysed by two different ap- 
proximations, (i) a mean-field approximation that neglects any local interaction in the network, 
and (ii) a local approximation considering correlations between pairs of nearest neighbors. In 
order to test the validity of these approximations, we compare their preditions with the averaged 
outcome of the microscopic computer simulations. We like to emphasize that our approach is 
general enough to be applied to various forms of frequency dependent processes on homogeneous 
networks with different number of neighbors. Even if a two-dimensional regular lattice is used to 
illustrate the dynamics, the approach is not restricted to that. 
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Our main result, in addition to the general framework for nonlinear frequency dependent pro- 
cesses, is the derivation of a phase diagram using the pair approximation derived in this paper. 
This approach predicts correctly both the type of the dynamics and the asymptotic values of the 
key variables dependent on the possible nonlinearities for the case of local interaction, m = 4. 
The predicted phase diagram was verified by comparision with extensive microscopic computer 
simulations rastering the whole parameter space. While the structure of the phase diagram was 
already known from previous computer simulations presented in [30J] we could demonstrate that 
the pair approximation works very well both for predicting the correct phase boundaries and the 
dynamics within these phases. It should be noticed that the pair approximation is a valuable 
tool, particularly with respect to computational efforts. The computer simulations are much more 
timeconsuming, since the results of the different runs have to be averaged afterwards. The pair 
approximation, on the other hand, is based on only 2 coupled equations and therefore needs less 
computational effort. In the following, we discuss some of the interesting findings. 

The region of correlated coexistence: Analysing the nonlinear VM with local interaction has 
shown that there are in fact only three different dynamic regimes dependent on the nonlinearities 
(ai,«2): (i) complete invasion, (ii) random coexistence, and (iii) correlated coexistence. The 
first one is already known as the standard behavior of the linear VM. Consequently the only 
interesting feature, namely the time to reach the ordered state dependent on the network size 
and topology, has been the subject of many investigations @, [48, 49]. Number (ii), on the other 
hand, only leads to trivial results as no real dynamics is observed. Thus, the most interesting 
regime is (iii) correlated coexistence, which can be found in a small, but not negligible parameter 
region below the voter point. This region separates the two dominant regimes (i) and (ii) and 
therefore was called a separatrix here. Going over from the right to the left side of the phase 
diagram in that region, we notice a transition from 0.5 to 1.0 (or 0.0 respectively) in the mean 
frequency, and from 0.5 to 0.7 to 1.0 in the local correlations. Thus, in fact c^i separates the two 
dynamic regimes (i) and (ii) (below the voter point). For parameters chosen from that region, we 
find in the CA simulations a long-term and nonstationary coexistence between the two species 
as on the right side of the phase diagram. But we also find the long-range spatial correlations 
that lead to the formation of spatial domains as shown e.g. in Fig. [6] - which is characteristic 
for the left side of the phase diagram. The spatial pattern formation is also indicated by large 
fluctuations of (x) shown in Fig. flTTtop). A single run, as shown in Fig. [7l clearly indicates the 
long-term nonstationary coexistence of both species. 

We emphasize that the separatrix between the two dynamic regimes is well predicted by the 
macroscopic dynamics resulting from the pair approximation (as can be clearly seen in Fig. [12]) . 
Most importantly, we could verify that the correlated coexistence of both species is not simply 
the effect of an additional perturbation, but results from the nonlinear interaction. 

Comparison with the mean-field phase diagram: In our paper, the mean-field approximation plays 
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the role of a reference state used to demonstrate the differences of the local analysis. The phase 
diagram of Fig. [3] distinguishes between six different regions, whereas the one in the local case, 
Fig. [12] (top) shows only three. Comparing the two phase diagrams, we realize that the most 
interesting regions in Fig. [31 namely (c) and (f), have simply collapsed into the separatrix shown 
in Fig. [T2]. The region (c) of unstable asymmetric coexistence or multiple outcome, respectively 
(see Sect. I4.2j l. relates to the separatrix line above the voter point. It should be noticed that the 
phase diagram for local interaction, Fig. [12], correctly predicts that the deterministic behavior 
for (01,02) = (0,1) leads to complete invasion (see Sect. 15.21 and Fig. 5(a)). 

The region (f) of stable asymmetric coexistence relates to the extended area of the separatrix 
shown below the voter point in Fig. [T21 where we still see a coexistence of both species - but the 
asymmetry between the two species relates to their changing dominance over time, as Figs.lTlfTTI 
(top) clearly illustrate. We conclude that in the local case no regions of stationary and asymmetric 
coexistence between the two species exist, as was predicted by the mean-field analysis. However, 
we find a (small but extended) region on the separatrix that shows the nonstationary and 
asymmetric coexistence of the two species for single realizations (which results in a symmetric 
coexistence averaged over runs, (x) = 0.5, see Fig. [TTJ top). 

The role of positive frequency dependence: The possible nonlinear responses in frequency depen- 
dent processes can be distinguished in four parameter areas of positive and negative frequency 
dependence and positive and negative allee effects, as Fig. [2] shows. Previous investigations [30J| 
assigned a dynamic leading to invasion to a positive frequency dependence, while associating a 
spatial coexistence with negative frequency dependence. Our investigations have shown that such 
an assignment does not unambiguously hold. In particular, a random coexistence can be found 
for negative frequency dependent dynamics as well as for positive frequency depencence, which 
was so far assigned to complete invasion only [3J][]. On the other hand, complete invasion is not 
observed only for positive frequency dependence, but also for positive and negative allee effects. 
A random spatial coexistence can be found for positive and negative allee dynamics as well. The 
only case where just one dynamic regime can be observed is the case of negative frequency de- 
pendence. We note, however, that the nonstationary long-term coexistence with spatial pattern 
formation occurs both for the positive frequency dependence and the negative allee dynamics, 
given that the parameters are chosen from the most interesting zone of the separatrix below the 
voter point. 

In conclusion, the region of positive frequency dependence bears in fact a much more richer 
dynamics, as it is transected by the separatrix we identified in the local analysis and thus shows 
all three types of dynamics we could identify for nonlinear voter models, namely (i) complete 
invasion of one of the species via the formation of large domains, (ii) long-term coexistence of 
both species with random distribution, (iii) long-term coexistence of both species with formation 
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of nonstationary domains. However, the most interesting regime (iii) is not restricted to positive 
frequency dependent processes, but can be also found for some negative allee effects. 

We summarize our findings by pointing out that nonlinear VM show indeed a very rich dynam- 
ics which was not much investigated yet. In addition to the phenomenon of complete invasion 
(or consensus) which occurs also beyond the linear VM, we find most interesting that certain 
parameter settings lead to a dynamics with nonstationarity and long-term correlations. Thinking 
about possible applications of the VM, we see that in particular this region has the potential to 
model relevant observations, be it the temporal dominance of certain species in a habitat or the 
temporal prevalence of certain opions (or political parties) in a social system. The nonstation- 
arity observed gives rise to the prediction that such dominance may not be the end, and change 
happens (even without additional perturbation). 
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Appendix 

Here we derive some explicit expressions for the three equations of the pair approximation dis- 
cussed in Sect. 13.31 for the global frequency (x) (eqn. (|3T1) ). the doublet frequency (x^i) and the 
correlation term cui, eqn. (l3~2|)). The equations are derived for the neighborhood m = 4. We use 
the notation x = (x). Using eqn. (|29l) and the transition rates of eqn. ([33]) . we find for (x), eqn. 



32/19] 



Frank Schweitzer, Laxmidhar Behera: 
Nonlinear Voter Models: The Transition from Invasion to Coexistence 

European Physical Journal B, vol 67 (2009). 



(|31l) in pair approximation: 
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We note that cm = c^g = x and c |o = c o|i = 1 — x in the mean-field limit, in which case eqn. 
(|38l) reduces to eqn. (|3ol) . 

In order to calculate the time derivative of the doublet frequency (xi t i) we have to consider how it 
is affected by changes of a in a specific occupation pattern of size m = 4, a° = {er, 02, 03, 04}, 
considering the <ij as constant. Again, in a frequency dependent process it is assumed that the 
transition does not depend on the exact distribution of the Oj, but only on the frequency of a 
particular state a in the neighborhood. Let S a> g describe a neighborhood where the center node 
in state a is surrounded by q nodes of the same state a. For any given q < m, there are (™) 
such occupation patterns. The global frequency of neighborhood S a>q is denoted as x a ^ q with the 
expectation value (x ajq ). Obviously, be calculated from the global frequencies x a j of 

all possible occupation distributions a' ( eqn. ©), that match the condition 



(39) 



i.e. it is defined as 



b cr,q 



E 



(40) 



Regarding the possible transitions, we are only interested in changes of the doublet (1,1), 
i.e. transitions (1,1) — ► (0,1) or (0,1) — > (1,1). The transition rates shall be denoted as 
io((0, 1)|(1, 1), S at q) and w((l, 1)|(0, 1), S a ,q) respectively, which of course depend on the local 
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neighborhood With this, the dynamics of the expected doublet frequency can be described 
by the rate equation: 

A m 

T (x 1>1 )(t)=J2[™[^ 1 WA),S , q ) (x , q ) 

9=0 (41) 

-w((0, 1)|(1, l),5i, 9 ) (x lig ) 

In order to specify the transition rates of the doublets to ((a 7 , 1) | (<r, 1), S a ,q), with a' = 1 — a 
and a 6 {0,1}, we note that there are only 10 distinct configurations of the neighborhood. Let 
us take the example a = {1,1,1,1,1}. A transition 1 — > of the center node would lead to 
the extinction of 4 doublets (a, (Tj) = (1, 1). On the other hand, the transition rate of the center 
node is e as known from eqn. (133"]) . This would result in to((0, 1)|(1, 1), 5*1,4) oc 4e. However, for 
a lattice of size N the number of doublets is 2N ', whereas there are exactly N neighborhoods a°. 
Therefore, if we apply the transition rates of the single nodes, eqn. (f33l) . to the transition of the 
doublets, their rates have to be scaled by 2. Similarly, if we take the example ct° = {0, 1, 1, 1,0}, 
a transition of the center node — > 1 would occur at the rate 1 — a<i and would create 3 new 
doublets. Applying the scaling factor of 2, we verify that w((l, 1)|(0, 1), 5o,i) = 3/2 (1 — c^)- 
This way we can determine the other possible transition rates: 
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Note that two of the transition rates are zero, because the respective doublets (1,1) or (0,1) do 

not exist in the assumed neighborhood. Finally, we express the (x CT)9 ) in eqn. (11TD by the < 

of eqn. (l4*0l) and apply the pair approximation, eqn. (|30l) . to the latter one. This way, we arrive 
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at the dynamic equation for (sci,i): 
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The third equation, eqn. (j32j) . for the correlation term ciw can be obtained in explicte form by 
using eqn. ([38]) for (x) and eqn. (l4*3l) for (2:1,1): 
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